11 January 2022
Read the taxonomy files back into Rstudio after performing the blastn search.
library(dplyr)
library(tidyr)
library(stringr)
library(ggplot2)
pool_taxonomy <-read.delim("blast_outputs/12S_poolseqs_taxlineage.txt", header = FALSE)
head(pool_taxonomy)
# use the full taxonomy rather than the seq id to collapse identical entries
pool_tax_df <- pool_taxonomy %>%
filter(V4 > 100) %>% # make sure all retained matches are >100 bp
select(-V2, -V5, -V6, -V7, -V14, -V8, -V9, -V10, -V11, -V12, -V15, -V16) %>% #remove unnecessary columns
group_by(V1, V17) %>% # group by the sequence key and the full taxonomy to reduce duplicate entries
unique() %>% # doing that reduced the number of entries from 146k to 17k
rename(qseqid=V1, perc_id=V3, length=V4, sscinames=V13, taxonomy=V17) %>% #rename headers
filter(!str_detect(sscinames, "environmental")) %>% # filter out any environmental samples
filter(!str_detect(sscinames, "synthetic")) %>% # filter out any synthetic "samples"
filter(perc_id >= 98) # seems like some of the matches below 98% are dubious (jellyfish and herring <1% different??)
Maybe break up that taxonomic info so that we can focus on the things that we can identify to some relevant level (at least order!?)
# first figure out where there are issues with that taxonomy format:
genus.hybrids <- pool_tax_df %>%
filter(!str_detect(taxonomy, ";")) %>%
# looking at these, they're probably hybrids or not good species, which is why the taxonomy broke down.
separate(sscinames, into = c("genus"), extra = "drop", remove = FALSE)
# and where there aren't issues
taxon_df <- pool_tax_df %>%
filter(str_detect(taxonomy, ";")) %>%
separate(taxonomy, into=c("kingdom", "phylum", "class", "order", "family", "genus", "species"), sep = ";")
# create a separate taxonomic lineage df with unique entries
lineage_df <- taxon_df %>%
ungroup() %>%
select(kingdom, phylum, class, order, family, genus) %>%
unique()
# there are only 18 unique genera in this small messed up part of the dataset
genus_fixed <- genus.hybrids %>%
ungroup() %>%
select(genus) %>%
unique() %>% # can I fix the lineage by adding existing data?
left_join(lineage_df, by = "genus") %>% # now the Tachyeres duck is the only loner
# fix that manually
mutate(family = ifelse(genus == "Tachyeres",
"Anatidae", family)) %>%
mutate(order = ifelse(genus == "Tachyeres",
"Anseriformes", order)) %>%
mutate(class = ifelse(genus == "Tachyeres",
"Aves", class)) %>%
mutate(phylum = ifelse(genus == "Tachyeres",
"Chordata", phylum)) %>%
mutate(kingdom = ifelse(genus == "Tachyeres",
"Eukaryota", kingdom))
# add the lineage info back to the rest of the data
genus_fixed_df <- genus.hybrids %>%
ungroup() %>%
left_join(genus_fixed, by = "genus") %>%
select(-taxonomy)
# and merge them together
fixed_tax_df <- taxon_df %>%
bind_rows(genus_fixed_df) %>%
# and fix the species syntax with underscores
mutate_if(is.character, str_replace_all, ' ', '_') %>%
filter(!str_detect(sscinames, "_x_"))
`mutate_if()` ignored the following grouping variables:
Column `qseqid`
# remove hybrids because the mtDNA sequence actually is only one or the other of the parent species
# how many total sequences?
fixed_tax_df %>%
select(qseqid) %>%
unique()
NA
There are four categories: 1. sequences that match a single species unambiguously (the minority)
Sequences that match multiple species are divided in three categories: 2. top matches > 2% identity than second-ranked matches 3. top matches < 2% identity than second-ranked matches 4. Multiple top matches with the same % identity
# 1. sequences that are unambiguously a single species
single_spp_seqs <- fixed_tax_df %>%
group_by(qseqid) %>%
add_tally(name = "n_taxa") %>%
filter(n_taxa == 1)
140 sequences are single-species
# remove the single-species seqs from the dataframe and then rank the hits by % identity for the remaining seqs
seq_id_diff <- fixed_tax_df %>%
anti_join(., single_spp_seqs) %>%
select(-length) %>%
group_by(qseqid, species, genus, family, order, class, phylum, kingdom) %>%
mutate(seq_percID = max(perc_id)) %>%
group_by(qseqid, species, genus, family, order, class, phylum, kingdom, seq_percID) %>%
summarise(max(seq_percID)) %>% # take just the top hit for each taxon (for each sequence)
select(-`max(seq_percID)`) %>%
ungroup() %>%
group_by(qseqid) %>%
mutate(id_rank = rank(desc(seq_percID), ties.method = "min")) %>% # rank the taxonomic hits per sequence by % id
mutate(top_perc = max(seq_percID)) %>% # designate the highest % id for the best taxonomic hit in each sequence (in some, but not all cases, this is 100%)
mutate(diff = top_perc - seq_percID) %>% # calculate the difference between the % identity of the top hit and each subsequent taxonomic hit
arrange(diff)
Joining, by = c("qseqid", "perc_id", "length", "sscinames", "kingdom", "phylum", "class", "order", "family", "genus", "species")
`summarise()` has grouped output by 'qseqid', 'species', 'genus', 'family', 'order', 'class', 'phylum', 'kingdom'. You can override using the `.groups` argument.
seq_id_diff %>%
filter(diff > 0)
NA
Now I have the single best entry for each species for each sequence ranked and with the difference between the first and second ranked entries calculated.
For sequences with multiple top hits, where the difference between ranked taxa = 0, I will end up defaulting to genus- or family-level ID (or carrying the individual species info around in some capacity). I will do the same for any sequences where the difference betweeen the first and second ranked taxa is < 2%.
Figure out which differences are > 2% and eliminate those first?
# filter out any taxa that are >2% less matching identity than the top taxonomic hit for a given sequence
to_remove_low_perc_hits <- seq_id_diff %>%
ungroup() %>%
group_by(qseqid) %>%
filter(diff > 2)
keepers <- seq_id_diff %>%
anti_join(to_remove_low_perc_hits)
Joining, by = c("qseqid", "species", "genus", "family", "order", "class", "phylum", "kingdom", "seq_percID", "id_rank", "top_perc", "diff")
# this data frame includes only those taxonomic hits that should be considered.
# so now I need to determine whether they should be assigned to genus, family, order, etc.
singletons <- keepers %>%
select(qseqid) %>%
tally() %>%
filter(n == 1)
# these are the seqs that now have only a single match
singleton_df <- singletons %>%
left_join(keepers) %>%
select(-n) %>%
bind_rows(single_spp_seqs) # combine the single spp data
## Genus-level matches
# remove the singletons from the bigger df
single_genus <- keepers %>%
anti_join(singleton_df)%>% # at best, these should be genus-level matches
group_by(qseqid, genus) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single genus
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "genus") %>%
mutate(taxon = genus)
## Family-level matches
single_family <- keepers %>%
anti_join(singleton_df)%>%
anti_join(single_genus) %>%
group_by(qseqid, family) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single family
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "family") %>%
mutate(taxon = family)
## Order-level matches
single_order <- keepers %>%
anti_join(singleton_df)%>%
anti_join(single_genus) %>%
anti_join(single_family) %>%
group_by(qseqid, order) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single order
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "order") %>%
mutate(taxon = order)
## Class-level matches
single_class <- keepers %>%
anti_join(singleton_df)%>%
anti_join(single_genus) %>%
anti_join(single_family) %>%
anti_join(single_order) %>%
group_by(qseqid, class) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single class
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "class") %>%
mutate(taxon = class)
## Phylum-level matches
single_phylum <- keepers %>%
anti_join(singleton_df)%>%
anti_join(single_genus) %>%
anti_join(single_family) %>%
anti_join(single_order) %>%
anti_join(single_class) %>%
group_by(qseqid, phylum) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single phylum
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "phylum") %>%
mutate(taxon = phylum)
# Kingdom level matches
single_kingdom <- keepers %>%
anti_join(singleton_df)%>%
anti_join(single_genus) %>%
anti_join(single_family) %>%
anti_join(single_order) %>%
anti_join(single_class) %>%
group_by(qseqid, kingdom) %>%
tally() %>%
ungroup() %>%
group_by(qseqid) %>%
tally() %>%
filter(n == 1) %>% # seqs that match a single kingdom
select(-n) %>%
left_join(., keepers) %>%
mutate(taxonomic_level = "kingdom") %>%
mutate(taxon = kingdom)
Seems… unlikely that a jellyfish would be 1% dissimilar from a herring. BUT, it’s below 98%, so maybe that’s at least one reason to be more stringent on the upstream side of this?
Adding that 98% filter in back up top.
Modify the singleton_df to include the right variable headers
single_spp <- singleton_df %>%
select(-perc_id, -length, -sscinames, -n_taxa) %>%
mutate(taxonomic_level = "species") %>%
mutate(taxon = species)
Okay, so that df should contain the basic taxa in the dataset assigned to the taxonomic level for which we can be relatively confident (>2% difference between the top-and-next highest hit).
I fundamentally have an issue with this step: we want to use biogeographic information to remove illogical matches, but we know that Fish Atlas is not comprehensive and so this “filter” is a very manual and tedious check. For example, sablefish, which is very common in nearshore SE AK, is not included in the Fish Atlas data.
Given that, here is how I will proceed:
Using family as the level of matching: First, I want to fill out the order for the fish atlas data using my data frame, then I’ll filter the data by family, but retain the species and genus designations as appropriate base on the prior taxonomic filtering step.
# read in fish atlas data from Pat:
fishatlas <- read.csv("../data/FAdata.csv") %>%
# remove unnecessary columns
select(Sp_CommonName, Sp_ScientificName, Fam_ScientificName) %>%
arrange(Fam_ScientificName)
# clean-up the format (remove spaces)
fishb <- fishatlas %>%
mutate_all(~(str_replace(., " ", "_"))) %>%
rename(common = Sp_CommonName, family = Fam_ScientificName)
# fix rockfishes so that it matches our eDNA data frame
fish.atlas.fixed <- fishb %>%
mutate(family = ifelse(str_detect(common, "rockfish"), "Sebastidae", family))
quick check: which taxa are filtered by the fish atlas?
# which taxa are filtered out when using the fish atlas?
# it turns out some of these are definitely real detections.
just.fish %>%
anti_join(., fish.atlas.fixed) %>%
unique() %>%
arrange(family) %>%
filter(!str_detect(species, "Carassius")) %>%
filter(!str_detect(species, "Cyprinus")) %>%
#select(family) %>%
unique() %>%
arrange(qseqid) %>%
filter(taxonomic_level %in% c("species", "genus", "family", "order")) %>%
select(qseqid, species, family, taxon, taxonomic_level)
Joining, by = "family"
**that is where sablefish and rockfish are filtered out! Why?? Manual checks on the taxa filtered out by Fish Atlas.
Focus on the “taxon” field here, since those are the data that will be carried forward The Merluccius should be filtered out because of geography. Lepomis_macrochirus = bluegill Labidesthes_sicculus = brook silverside Leuciscidae = freshwater fish (Cyprinidae) Ambloplites = freshwater fish Merluccius = Atlantic/Europe Mallotus_villosus = capelin (possible in AK??) Micropterus_dolomieu = smallmouth bass Etheostoma = freshwater Engraulidae = most likely from the fish food Lipolagus_ochotensis = black-eared smelt (deep water… but in the N Pacific) Catostomidae = freshwater suckers Percidae = fresh/brackish waters of N America Anoplopomatidae = sablefish!!! Scomber_scombrus = Atlantic mackerel (fish food, maybe) Echiichthys_vipera = Atlantic Micropterus = black basses Macrouridae = deep water (Grenadiers or rattails) Pimephales = freshwater Percina_peltata = freshwater/Eastern US Trachinus_draco = Atlantic Bagre_marinus = catfish Anarrhichthys_ocellatus = wolf eel
Here’s where the painful part comes in - we want to keep the following species that are otherwise filtered out by the Fish Atlas join:
“Anoplopomatidae” “Anarrhichthys_ocellatus” “Macrouridae” “Lipolagus_ochotensis” “Mallotus_villosus”
Probably the simplest way to do this is to keep the inner join with Fish Atlas and then add these back on.
fishatlas.edna.tax <- fish.atlas.fixed %>%
inner_join(., just.fish, by = "family") %>%
mutate(to.filter = ifelse(taxonomic_level == "species" & species != Sp_ScientificName, "yes", "no")) # this provides us a way to filter mis-identified taxa to the species level based on fish atlas data - but these taxa are fine to retain at the genus/family level
Quick question: how many taxa get removed by because the family matches but the species does not?
# remove duplicates so that I can take a closer look at which species-level IDs would be removed.
fishatlas.edna.tax %>%
filter(to.filter == "yes") %>%
select(-qseqid) %>%
select(taxon) %>%
unique()
NA
NA
In total, 8 species would be removed from the dataset due to matching family but mismatching species with Fish Atlas. Let’s verify each one.
Oncorhynchus_tshawytscha - Chinook salmon… definitely present in the marine nearshore outside of Juneau. Probably shouldn’t remove. Leptocottus_armatus - Pacific staghorn sculpin - also quite likely to be present Artedius_fenestralis - Padded sculpin - distribution extends to AK.. unknownn. Atheresthes_evermanni - Kamchatka flounder - unlikely Sardina_pilchardus - European pilchard - unlikely (fish food) Opisthonema_oglinum - Atlantic thread herring - unlikely (fish food) Hexagrammos_decagrammus - kelp greenling (definitely possible) Gymnammodytes_semisquamatus - Smooth sandeel - distribution = Europe (maybe fish food)
So the conclusion is that 50% of those are possible or likely in our study system and we should keep the following: Oncorhynchus_tshawytscha Leptocottus_armatus Artedius_fenestralis Hexagrammos_decagrammus Atheresthes_evermanni -should be switched to genus-level because arrowtooth flounder is a sister-species of Kamchatka flounder
I hate doing this, but I’m going to have to manually filter the four to preserve the others:
# filter based on fish atlas, but retain species that were double-checked
fishatlas.edna.filtered <- fishatlas.edna.tax %>%
filter(!taxon %in% c("Sardina_pilchardus", "Opisthonema_oglinum", "Gymnammodytes_semisquamatus")) %>%
# and then add back on the taxa from extra.keepers
bind_rows(extra.keepers) %>%
# switch Atheresthes_evermanni to a genus-level match
mutate(taxon = ifelse(species == "Atheresthes_evermanni",
"Atheresthes", taxon)) %>%
mutate(taxonomic_level = ifelse(species == "Atheresthes_evermanni",
"genus", taxonomic_level)) %>%
select(-to.filter)
And just to confirm that that took care of the outstanding issues, look for sebastes, sablefish, etc.
# export that fish-only taxonomy df
fishatlas.edna.filtered %>%
write.csv("csv_outputs/12Sfishatlas.edna.taxonomy.csv")
There are only 748 unique sequences that made it through the blast/taxonomy gauntlet, so why are there an order of magnitude more sequences in the seqtab.nochim data tables. Maybe because they didn’t make it through the BLAST filtering? i.e., no match at >98% identity?
go to 06-merge-samples-w-taxonomy.Rmd instead for the latest version of this analysis.
asv_tbl <- read.csv("csv_outputs/12S_poolseq_ASVtable.csv") %>%
rename(sample.PCR = X) %>%
separate(sample.PCR, into = c("sample", "rep"), remove = F)
Combine taxonomy and sample ASV table
# pivot longer for ease of joining
asv_long_tbl <- asv_tbl %>%
pivot_longer(cols = 4:7705, names_to = "ASV", values_to = "count")
# join data frames
uncollapsed_poolseq_tax_df <- asv_long_tbl %>%
left_join(., sorted_tax_df, by = c("ASV" = "qseqid"))
# collapse these based on the appropriate taxonomic level
collapsed_poolseq_tax_df <- uncollapsed_poolseq_tax_df %>%
select(sample.PCR, sample, rep, ASV, count, taxon, taxonomic_level) %>% ## I need to go back and fix the missing % for the species-level hits
unique()
That data frame still contains a bunch of NAs that didn’t have blast hits. Let’s go ahead and remove them.
poolseq_tax_tidydf <- collapsed_poolseq_tax_df %>%
filter(!is.na(taxon)) %>% # remove NAs
filter(count > 0) # remove entries with zero reads for that sequence
Add metadata to that tidy frame
# metadata tidied and ready to go
meta <- read.csv("../metadata/amalga_clean_metadata.csv")
poolseq_tax_tidymeta <- poolseq_tax_tidydf %>%
left_join(., meta, by = "sample") #%>%
#filter(!is.na(label)) # remove aquarium samples?
poolseq_tax_filtered_w_meta <- poolseq_tax_tidymeta %>%
arrange(distance) %>%
filter(taxonomic_level != "kingdom") %>% # remove these totally uninformative matches
# collapse multiple entries for the same taxon in the same sample replicate
group_by(sample.PCR, taxon) %>%
mutate(count = sum(count)) %>%
select(-ASV) %>%
unique() # that reduces the number of rows from 7k to 4k
That’s pretty useful.
Let’s break it down in a few different ways:
# negatives?
poolseq_tax_filtered_w_meta %>%
filter(label == "blank", taxon != "Homo_sapiens") %>%
ggplot(aes(x=sample.PCR, y=count, fill=taxon)) +
geom_bar(stat="identity") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90),
axis.title.x = element_text(margin = margin(t=10))
) +
labs(
y = "reads",
x = "extraction blanks"
)
ggsave("pdf_outputs/12S_extraction_blanks.pdf")
Saving 7.26 x 4.49 in image
# pcr negatives
poolseq_tax_filtered_w_meta %>%
filter(str_detect(sample.PCR, "NEGATIVE"), taxon != "Homo_sapiens") %>%
ggplot(aes(x=sample.PCR, y=count, fill=taxon)) +
geom_bar(stat="identity") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90),
axis.title.x = element_text(margin = margin(t=10))
) +
labs(
y = "reads",
x = "PCR negative controls"
)
ggsave("pdf_outputs/pcr_negatives.pdf")
Saving 7.33 x 4.52 in image
# these PCR negatives look awesome.
# positives?
poolseq_tax_filtered_w_meta %>%
filter(str_detect(sample.PCR, "POSITIVE"), taxon != "Homo_sapiens") %>%
ggplot(aes(x=sample.PCR, y=count, fill=taxon)) +
geom_bar(stat="identity") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90),
axis.title.x = element_text(margin = margin(t=10))
) +
labs(
y = "reads",
x = "PCR positive controls"
)
I guess the goldfish doesn’t show up at any taxonomic level higher than class?
# species/genus/family-level matches?
poolseq_tax_filtered_w_meta %>%
filter(!is.na(label),
taxon != "Homo_sapiens",
taxon != "Canis_lupus",
taxonomic_level %in% c("species")) %>%
ggplot(aes(x=sample, y=count, fill=taxon)) +
geom_bar(stat="identity") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90),
axis.title.x = element_text(margin = margin(t=10))
) +
labs(
y = "reads",
x = "sample replicates"
)
NA
NA
Next up: organize by distance and look at composition. Also - could summarize/further collapse taxa by sample when sequences fall in multiple taxonomic levels.
library(vegan)
Loading required package: permute
Loading required package: lattice
This is vegan 2.5-7
To run an ordination you will need a data-frame consisting of plot by species (or trait) matrix AND a “groups” data-frame which should consist of plots with a coding variable for what group each plot belongs to, this will be used for plotting the ordination.
This is not really the format that we currently have, so let’s move things around a bit.
# we don't need the metadata for this - simplify
pool.vegan.df <- poolseq_tax_filtered_w_meta %>%
filter(!taxon %in% c("Homo_sapiens", "Canis_lupus")) %>% # remove these because
group_by(location, label, type) %>%
mutate(grp = group_indices()) %>% # this creates a group index for each set of field/pcr replicates
ungroup() %>%
select(grp, sample.PCR, taxon, count)
# separate out the grp info
grps <- pool.vegan.df %>%
select(grp, sample.PCR)
# separate out the df to convert
wide_df <- pool.vegan.df %>%
select(-grp) %>%
pivot_wider(id_cols=sample.PCR, names_from=taxon, values_from=count) %>%
replace(is.na(.), 0) %>% # replace NAs with zeros
select(-sample.PCR)
matrix <- data.matrix(wide_df, rownames.force = T)
convert absolute abundance to relative abundance
# Running NMDS in vegan (metaMDS)
NMS <-
metaMDS(tax.distmat,
distance = "bray",
k = 3,
maxit = 999,
trymax = 500,
wascores = TRUE)
Run 0 stress 0.1047505
Run 1 stress 0.1091038
Run 2 stress 0.1033705
... New best solution
... Procrustes: rmse 0.01517027 max resid 0.1431347
Run 3 stress 0.1031345
... New best solution
... Procrustes: rmse 0.01476544 max resid 0.169073
Run 4 stress 0.1099409
Run 5 stress 0.1075893
Run 6 stress 0.1034186
... Procrustes: rmse 0.0145414 max resid 0.1562133
Run 7 stress 0.1030296
... New best solution
... Procrustes: rmse 0.01489652 max resid 0.1528805
Run 8 stress 0.1035558
Run 9 stress 0.1039788
Run 10 stress 0.1028844
... New best solution
... Procrustes: rmse 0.01462787 max resid 0.1633978
Run 11 stress 0.105363
Run 12 stress 0.1054943
Run 13 stress 0.1111475
Run 14 stress 0.1053857
Run 15 stress 0.1058309
Run 16 stress 0.1029229
... Procrustes: rmse 0.01399598 max resid 0.1420114
Run 17 stress 0.110239
Run 18 stress 0.1036202
Run 19 stress 0.1050009
Run 20 stress 0.103275
... Procrustes: rmse 0.01339265 max resid 0.164441
Run 21 stress 0.103434
Run 22 stress 0.1040243
Run 23 stress 0.1044319
Run 24 stress 0.1122533
Run 25 stress 0.109398
Run 26 stress 0.1073117
Run 27 stress 0.1053562
Run 28 stress 0.1107769
Run 29 stress 0.1043899
Run 30 stress 0.1112774
Run 31 stress 0.1095436
Run 32 stress 0.1038537
Run 33 stress 0.1115517
Run 34 stress 0.1054331
Run 35 stress 0.1136744
Run 36 stress 0.1071017
Run 37 stress 0.106483
Run 38 stress 0.1033172
... Procrustes: rmse 0.01208219 max resid 0.1029581
Run 39 stress 0.1040138
Run 40 stress 0.1031196
... Procrustes: rmse 0.01269745 max resid 0.1236347
Run 41 stress 0.1037689
Run 42 stress 0.1040558
Run 43 stress 0.1023241
... New best solution
... Procrustes: rmse 0.0108321 max resid 0.1406588
Run 44 stress 0.1075474
Run 45 stress 0.1042029
Run 46 stress 0.1075343
Run 47 stress 0.1036712
Run 48 stress 0.1044605
Run 49 stress 0.1093957
Run 50 stress 0.1047364
Run 51 stress 0.1071635
Run 52 stress 0.1040558
Run 53 stress 0.1069962
Run 54 stress 0.1078106
Run 55 stress 0.1035598
Run 56 stress 0.1032626
Run 57 stress 0.1042423
Run 58 stress 0.1104173
Run 59 stress 0.1068436
Run 60 stress 0.1037441
Run 61 stress 0.1066773
Run 62 stress 0.1064878
Run 63 stress 0.1139548
Run 64 stress 0.108836
Run 65 stress 0.1031774
Run 66 stress 0.1035821
Run 67 stress 0.1041209
Run 68 stress 0.1110943
Run 69 stress 0.1066009
Run 70 stress 0.1110221
Run 71 stress 0.1089743
Run 72 stress 0.1136639
Run 73 stress 0.107876
Run 74 stress 0.1046272
Run 75 stress 0.1141784
Run 76 stress 0.1105466
Run 77 stress 0.1040599
Run 78 stress 0.1131514
Run 79 stress 0.1081336
Run 80 stress 0.1058468
Run 81 stress 0.1072477
Run 82 stress 0.1032487
Run 83 stress 0.1124266
Run 84 stress 0.1053103
Run 85 stress 0.1128828
Run 86 stress 0.1033269
Run 87 stress 0.1073788
Run 88 stress 0.1075338
Run 89 stress 0.1067032
Run 90 stress 0.1033845
Run 91 stress 0.10491
Run 92 stress 0.1053343
Run 93 stress 0.1038026
Run 94 stress 0.1039996
Run 95 stress 0.1057072
Run 96 stress 0.1043389
Run 97 stress 0.1031114
Run 98 stress 0.1085936
Run 99 stress 0.1035071
Run 100 stress 0.1085043
Run 101 stress 0.1043153
Run 102 stress 0.104713
Run 103 stress 0.1077827
Run 104 stress 0.1062981
Run 105 stress 0.1031222
Run 106 stress 0.1084776
Run 107 stress 0.1057907
Run 108 stress 0.1039365
Run 109 stress 0.1091723
Run 110 stress 0.1080755
Run 111 stress 0.1034593
Run 112 stress 0.1061215
Run 113 stress 0.1108576
Run 114 stress 0.1059668
Run 115 stress 0.1055788
Run 116 stress 0.1115548
Run 117 stress 0.1043042
Run 118 stress 0.1042351
Run 119 stress 0.1074185
Run 120 stress 0.103832
Run 121 stress 0.1032032
Run 122 stress 0.1099107
Run 123 stress 0.1068055
Run 124 stress 0.1045153
Run 125 stress 0.1058854
Run 126 stress 0.1085682
Run 127 stress 0.1034648
Run 128 stress 0.1053297
Run 129 stress 0.1082581
Run 130 stress 0.103385
Run 131 stress 0.1038871
Run 132 stress 0.1031488
Run 133 stress 0.1099827
Run 134 stress 0.1115565
Run 135 stress 0.10442
Run 136 stress 0.1040222
Run 137 stress 0.10662
Run 138 stress 0.1056669
Run 139 stress 0.1067992
Run 140 stress 0.106493
Run 141 stress 0.1083335
Run 142 stress 0.1067389
Run 143 stress 0.1070057
Run 144 stress 0.1079087
Run 145 stress 0.1059276
Run 146 stress 0.103159
Run 147 stress 0.1092317
Run 148 stress 0.1073626
Run 149 stress 0.1095654
Run 150 stress 0.1068984
Run 151 stress 0.1079596
Run 152 stress 0.1110705
Run 153 stress 0.1041155
Run 154 stress 0.1042449
Run 155 stress 0.1137419
Run 156 stress 0.1120904
Run 157 stress 0.1068519
Run 158 stress 0.1026228
... Procrustes: rmse 0.008829727 max resid 0.1013661
Run 159 stress 0.1100415
Run 160 stress 0.1044317
Run 161 stress 0.10397
Run 162 stress 0.1079679
Run 163 stress 0.105392
Run 164 stress 0.1040166
Run 165 stress 0.1043113
Run 166 stress 0.1038135
Run 167 stress 0.1030774
Run 168 stress 0.1028823
Run 169 stress 0.1104579
Run 170 stress 0.1105455
Run 171 stress 0.1032114
Run 172 stress 0.1088956
Run 173 stress 0.1080987
Run 174 stress 0.1042083
Run 175 stress 0.1040545
Run 176 stress 0.1055706
Run 177 stress 0.1078103
Run 178 stress 0.1039281
Run 179 stress 0.104422
Run 180 stress 0.1038253
Run 181 stress 0.1050384
Run 182 stress 0.1060571
Run 183 stress 0.1089949
Run 184 stress 0.1052862
Run 185 stress 0.1040233
Run 186 stress 0.1131226
Run 187 stress 0.1100965
Run 188 stress 0.1065048
Run 189 stress 0.1041713
Run 190 stress 0.1123618
Run 191 stress 0.1041876
Run 192 stress 0.1067663
Run 193 stress 0.1029046
Run 194 stress 0.1080911
Run 195 stress 0.1086436
Run 196 stress 0.1109009
Run 197 stress 0.1078054
Run 198 stress 0.1092383
Run 199 stress 0.113815
Run 200 stress 0.1102316
Run 201 stress 0.1053111
Run 202 stress 0.1077149
Run 203 stress 0.1078718
Run 204 stress 0.1041129
Run 205 stress 0.1097915
Run 206 stress 0.1066643
Run 207 stress 0.1035781
Run 208 stress 0.1045359
Run 209 stress 0.1087377
Run 210 stress 0.1071641
Run 211 stress 0.1132297
Run 212 stress 0.1130048
Run 213 stress 0.1056566
Run 214 stress 0.1122072
Run 215 stress 0.109456
Run 216 stress 0.1037874
Run 217 stress 0.1073794
Run 218 stress 0.1035854
Run 219 stress 0.1061877
Run 220 stress 0.108798
Run 221 stress 0.1062357
Run 222 stress 0.1040098
Run 223 stress 0.1132609
Run 224 stress 0.114126
Run 225 stress 0.1065758
Run 226 stress 0.1059224
Run 227 stress 0.1041753
Run 228 stress 0.1033392
Run 229 stress 0.1034376
Run 230 stress 0.1112893
Run 231 stress 0.1034152
Run 232 stress 0.1147076
Run 233 stress 0.1071721
Run 234 stress 0.1042944
Run 235 stress 0.1050444
Run 236 stress 0.108261
Run 237 stress 0.1065422
Run 238 stress 0.1023964
... Procrustes: rmse 0.01074999 max resid 0.1310516
Run 239 stress 0.1043985
Run 240 stress 0.108586
Run 241 stress 0.1065348
Run 242 stress 0.1071631
Run 243 stress 0.110684
Run 244 stress 0.1083524
Run 245 stress 0.1053337
Run 246 stress 0.1078528
Run 247 stress 0.104499
Run 248 stress 0.103487
Run 249 stress 0.1047615
Run 250 stress 0.1035198
Run 251 stress 0.103747
Run 252 stress 0.1058733
Run 253 stress 0.1034251
Run 254 stress 0.1067844
Run 255 stress 0.1086243
Run 256 stress 0.1045882
Run 257 stress 0.1063707
Run 258 stress 0.1030918
Run 259 stress 0.1040189
Run 260 stress 0.1082816
Run 261 stress 0.1095282
Run 262 stress 0.1059059
Run 263 stress 0.1084944
Run 264 stress 0.1046809
Run 265 stress 0.1096994
Run 266 stress 0.1093923
Run 267 stress 0.1028253
Run 268 stress 0.1149591
Run 269 stress 0.1047812
Run 270 stress 0.1078624
Run 271 stress 0.105012
Run 272 stress 0.1125052
Run 273 stress 0.1057068
Run 274 stress 0.1092531
Run 275 stress 0.104797
Run 276 stress 0.1031935
Run 277 stress 0.1078094
Run 278 stress 0.1043631
Run 279 stress 0.1092091
Run 280 stress 0.1044983
Run 281 stress 0.1041113
Run 282 stress 0.1053773
Run 283 stress 0.1034588
Run 284 stress 0.1083526
Run 285 stress 0.1055744
Run 286 stress 0.1036517
Run 287 stress 0.1043398
Run 288 stress 0.1053034
Run 289 stress 0.1067655
Run 290 stress 0.105802
Run 291 stress 0.107568
Run 292 stress 0.1103305
Run 293 stress 0.112835
Run 294 stress 0.1053336
Run 295 stress 0.1037726
Run 296 stress 0.1029617
Run 297 stress 0.1033539
Run 298 stress 0.1059967
Run 299 stress 0.1110852
Run 300 stress 0.1046125
Run 301 stress 0.1064393
Run 302 stress 0.1080683
Run 303 stress 0.1063069
Run 304 stress 0.1054926
Run 305 stress 0.1116901
Run 306 stress 0.103548
Run 307 stress 0.1081589
Run 308 stress 0.1125721
Run 309 stress 0.1097499
Run 310 stress 0.1039662
Run 311 stress 0.1147889
Run 312 stress 0.1046511
Run 313 stress 0.1049009
Run 314 stress 0.1086474
Run 315 stress 0.1072537
Run 316 stress 0.1032225
Run 317 stress 0.1032559
Run 318 stress 0.1037649
Run 319 stress 0.1054288
Run 320 stress 0.1104593
Run 321 stress 0.107653
Run 322 stress 0.110044
Run 323 stress 0.1051683
Run 324 stress 0.1029237
Run 325 stress 0.1060015
Run 326 stress 0.1123066
Run 327 stress 0.1143683
Run 328 stress 0.1078046
Run 329 stress 0.1082902
Run 330 stress 0.1055701
Run 331 stress 0.1128816
Run 332 stress 0.1035831
Run 333 stress 0.1087047
Run 334 stress 0.1041095
Run 335 stress 0.1034921
Run 336 stress 0.1131037
Run 337 stress 0.1052898
Run 338 stress 0.1090951
Run 339 stress 0.1095752
Run 340 stress 0.1091098
Run 341 stress 0.1047604
Run 342 stress 0.1090653
Run 343 stress 0.1047757
Run 344 stress 0.1052119
Run 345 stress 0.1056775
Run 346 stress 0.1064509
Run 347 stress 0.1060602
Run 348 stress 0.1031868
Run 349 stress 0.1053675
Run 350 stress 0.1064393
Run 351 stress 0.1037076
Run 352 stress 0.1047263
Run 353 stress 0.1077758
Run 354 stress 0.1028465
Run 355 stress 0.1039849
Run 356 stress 0.1124709
Run 357 stress 0.1051077
Run 358 stress 0.1029264
Run 359 stress 0.1031474
Run 360 stress 0.1056758
Run 361 stress 0.1032205
Run 362 stress 0.1089219
Run 363 stress 0.1055382
Run 364 stress 0.1030962
Run 365 stress 0.1029831
Run 366 stress 0.1110793
Run 367 stress 0.1091858
Run 368 stress 0.1037874
Run 369 stress 0.1032391
Run 370 stress 0.1078163
Run 371 stress 0.1069993
Run 372 stress 0.1116935
Run 373 stress 0.1043484
Run 374 stress 0.1043161
Run 375 stress 0.1086366
Run 376 stress 0.1041711
Run 377 stress 0.1046042
Run 378 stress 0.1097481
Run 379 stress 0.1080457
Run 380 stress 0.1043261
Run 381 stress 0.103593
Run 382 stress 0.1047723
Run 383 stress 0.1050301
Run 384 stress 0.1081724
Run 385 stress 0.1137157
Run 386 stress 0.1048735
Run 387 stress 0.1060884
Run 388 stress 0.1059132
Run 389 stress 0.103554
Run 390 stress 0.1102933
Run 391 stress 0.1060102
Run 392 stress 0.1048778
Run 393 stress 0.108088
Run 394 stress 0.1124055
Run 395 stress 0.1073617
Run 396 stress 0.1131979
Run 397 stress 0.1058047
Run 398 stress 0.1035497
Run 399 stress 0.1050533
Run 400 stress 0.1060043
Run 401 stress 0.1032813
Run 402 stress 0.103049
Run 403 stress 0.1047513
Run 404 stress 0.1063919
Run 405 stress 0.1062946
Run 406 stress 0.1099276
Run 407 stress 0.1050615
Run 408 stress 0.1040607
Run 409 stress 0.1038954
Run 410 stress 0.1067165
Run 411 stress 0.1093918
Run 412 stress 0.1041121
Run 413 stress 0.1063597
Run 414 stress 0.1060783
Run 415 stress 0.1045171
Run 416 stress 0.1041195
Run 417 stress 0.1105354
Run 418 stress 0.1082509
Run 419 stress 0.1095471
Run 420 stress 0.1064907
Run 421 stress 0.106217
Run 422 stress 0.1044849
Run 423 stress 0.1044891
Run 424 stress 0.1039686
Run 425 stress 0.1034413
Run 426 stress 0.1065345
Run 427 stress 0.1051907
Run 428 stress 0.105302
Run 429 stress 0.103088
Run 430 stress 0.1078367
Run 431 stress 0.1081468
Run 432 stress 0.103194
Run 433 stress 0.1085562
Run 434 stress 0.1101137
Run 435 stress 0.1076278
Run 436 stress 0.1041889
Run 437 stress 0.1036967
Run 438 stress 0.1035539
Run 439 stress 0.1097133
Run 440 stress 0.1060898
Run 441 stress 0.1042543
Run 442 stress 0.1085263
Run 443 stress 0.1035459
Run 444 stress 0.1078768
Run 445 stress 0.1079541
Run 446 stress 0.1031527
Run 447 stress 0.1126149
Run 448 stress 0.103595
Run 449 stress 0.1112164
Run 450 stress 0.106703
Run 451 stress 0.1051164
Run 452 stress 0.1058897
Run 453 stress 0.1037466
Run 454 stress 0.1048177
Run 455 stress 0.1124754
Run 456 stress 0.1081922
Run 457 stress 0.1066093
Run 458 stress 0.1041059
Run 459 stress 0.1039044
Run 460 stress 0.1035431
Run 461 stress 0.1072625
Run 462 stress 0.1080013
Run 463 stress 0.1083675
Run 464 stress 0.1041363
Run 465 stress 0.1043774
Run 466 stress 0.1028718
Run 467 stress 0.1057578
Run 468 stress 0.1043604
Run 469 stress 0.1041471
Run 470 stress 0.1079567
Run 471 stress 0.1085731
Run 472 stress 0.1084608
Run 473 stress 0.1071668
Run 474 stress 0.1082753
Run 475 stress 0.1040823
Run 476 stress 0.1052025
Run 477 stress 0.1109382
Run 478 stress 0.1045455
Run 479 stress 0.1064157
Run 480 stress 0.1027789
... Procrustes: rmse 0.01328367 max resid 0.1524937
Run 481 stress 0.1081466
Run 482 stress 0.1088787
Run 483 stress 0.1037462
Run 484 stress 0.1078385
Run 485 stress 0.1027742
... Procrustes: rmse 0.01174669 max resid 0.1269963
Run 486 stress 0.1043833
Run 487 stress 0.108456
Run 488 stress 0.1033895
Run 489 stress 0.1139604
Run 490 stress 0.1080004
Run 491 stress 0.1064829
Run 492 stress 0.1037161
Run 493 stress 0.1039652
Run 494 stress 0.1041698
Run 495 stress 0.1099811
Run 496 stress 0.1081474
Run 497 stress 0.1047244
Run 498 stress 0.1119344
Run 499 stress 0.1025052
... Procrustes: rmse 0.01043141 max resid 0.09655726
Run 500 stress 0.1038412
*** No convergence -- monoMDS stopping criteria:
176: stress ratio > sratmax
324: scale factor of the gradient < sfgrmin
It is common for NMDS analyses to start by running with 2-dimensions (k), but you want to increase the number of dimensions to ensure a minimized stress value. Keep in mind that anything more than 5-dimensions makes it difficult to interpret a 2-dimensional plot.
As a rule of thumb literature has identified the following cut-off values for stress-level:
Higher than 0.2 is poor (risks for false interpretation). 0.1 - 0.2 is fair (some distances can be misleading for interpretation). 0.05 - 0.1 is good (can be confident in inferences from plot). Less than 0.05 is excellent (this can be rare).
A few options: 1) increase k 2) increase maxit
let’s try increasing k
# Running NMDS in vegan (metaMDS)
NMS.k4 <-
metaMDS(tax.distmat,
distance = "bray",
k = 4,
maxit = 999,
trymax = 500,
wascores = TRUE)
Run 0 stress 0.07386974
Run 1 stress 0.07491524
Run 2 stress 0.07437815
Run 3 stress 0.07403766
... Procrustes: rmse 0.01514743 max resid 0.1544086
Run 4 stress 0.07476686
Run 5 stress 0.07504459
Run 6 stress 0.07497444
Run 7 stress 0.0740845
... Procrustes: rmse 0.01147333 max resid 0.1218068
Run 8 stress 0.07525555
Run 9 stress 0.07520172
Run 10 stress 0.07445861
Run 11 stress 0.07519227
Run 12 stress 0.07405646
... Procrustes: rmse 0.01127032 max resid 0.128445
Run 13 stress 0.07502137
Run 14 stress 0.07589977
Run 15 stress 0.07495102
Run 16 stress 0.07395249
... Procrustes: rmse 0.01506798 max resid 0.1277471
Run 17 stress 0.07414531
... Procrustes: rmse 0.01108162 max resid 0.1278024
Run 18 stress 0.07411322
... Procrustes: rmse 0.01049023 max resid 0.1345969
Run 19 stress 0.07408317
... Procrustes: rmse 0.01416774 max resid 0.1311196
Run 20 stress 0.07425687
... Procrustes: rmse 0.008242095 max resid 0.1099863
Run 21 stress 0.07498896
Run 22 stress 0.07497087
Run 23 stress 0.07466173
Run 24 stress 0.0744282
Run 25 stress 0.07432908
... Procrustes: rmse 0.01199421 max resid 0.1274044
Run 26 stress 0.07414727
... Procrustes: rmse 0.0121427 max resid 0.1257852
Run 27 stress 0.0738444
... New best solution
... Procrustes: rmse 0.01148858 max resid 0.1131865
Run 28 stress 0.07412599
... Procrustes: rmse 0.01206593 max resid 0.127347
Run 29 stress 0.07386381
... Procrustes: rmse 0.01037911 max resid 0.1278972
Run 30 stress 0.07464
Run 31 stress 0.07401148
... Procrustes: rmse 0.01283616 max resid 0.1251787
Run 32 stress 0.07493285
Run 33 stress 0.07453871
Run 34 stress 0.07408212
... Procrustes: rmse 0.00942815 max resid 0.1120016
Run 35 stress 0.07491675
Run 36 stress 0.07408543
... Procrustes: rmse 0.01112055 max resid 0.1149321
Run 37 stress 0.0753959
Run 38 stress 0.0750333
Run 39 stress 0.07411234
... Procrustes: rmse 0.0122264 max resid 0.1155168
Run 40 stress 0.07430884
... Procrustes: rmse 0.01533867 max resid 0.1551778
Run 41 stress 0.07526799
Run 42 stress 0.07397533
... Procrustes: rmse 0.01238658 max resid 0.1158422
Run 43 stress 0.07423644
... Procrustes: rmse 0.01243427 max resid 0.1256216
Run 44 stress 0.07462486
Run 45 stress 0.07535054
Run 46 stress 0.07441149
Run 47 stress 0.07447003
Run 48 stress 0.07415679
... Procrustes: rmse 0.01289608 max resid 0.1260629
Run 49 stress 0.07449526
Run 50 stress 0.07419937
... Procrustes: rmse 0.01292097 max resid 0.1286732
Run 51 stress 0.07500998
Run 52 stress 0.07336348
... New best solution
... Procrustes: rmse 0.01031954 max resid 0.1285916
Run 53 stress 0.0742318
Run 54 stress 0.07603917
Run 55 stress 0.07387452
Run 56 stress 0.07489372
Run 57 stress 0.07403113
Run 58 stress 0.07486089
Run 59 stress 0.07417056
Run 60 stress 0.07422359
Run 61 stress 0.07418035
Run 62 stress 0.07410374
Run 63 stress 0.07430432
Run 64 stress 0.07485444
Run 65 stress 0.07444278
Run 66 stress 0.07511817
Run 67 stress 0.07379722
... Procrustes: rmse 0.01075407 max resid 0.1286048
Run 68 stress 0.07528454
Run 69 stress 0.07419061
Run 70 stress 0.07373564
... Procrustes: rmse 0.01101958 max resid 0.1269649
Run 71 stress 0.07502309
Run 72 stress 0.07490364
Run 73 stress 0.07423336
Run 74 stress 0.07413461
Run 75 stress 0.07469098
Run 76 stress 0.0751561
Run 77 stress 0.07394234
Run 78 stress 0.0751409
Run 79 stress 0.07397217
Run 80 stress 0.07350773
... Procrustes: rmse 0.008049176 max resid 0.1321206
Run 81 stress 0.07534628
Run 82 stress 0.07539287
Run 83 stress 0.07513191
Run 84 stress 0.07465574
Run 85 stress 0.07536908
Run 86 stress 0.07549339
Run 87 stress 0.07519677
Run 88 stress 0.07419013
Run 89 stress 0.07423067
Run 90 stress 0.07440885
Run 91 stress 0.07479592
Run 92 stress 0.07367861
... Procrustes: rmse 0.008423985 max resid 0.1316586
Run 93 stress 0.07577634
Run 94 stress 0.073537
... Procrustes: rmse 0.006278642 max resid 0.09073157
Run 95 stress 0.07456939
Run 96 stress 0.0742466
Run 97 stress 0.07426099
Run 98 stress 0.07533133
Run 99 stress 0.07511749
Run 100 stress 0.07384193
... Procrustes: rmse 0.01201916 max resid 0.1300418
Run 101 stress 0.0750367
Run 102 stress 0.0746003
Run 103 stress 0.07443359
Run 104 stress 0.07506973
Run 105 stress 0.07496083
Run 106 stress 0.0744362
Run 107 stress 0.07539374
Run 108 stress 0.07425111
Run 109 stress 0.07394821
Run 110 stress 0.07538046
Run 111 stress 0.07498602
Run 112 stress 0.07371408
... Procrustes: rmse 0.009580648 max resid 0.1308994
Run 113 stress 0.07369853
... Procrustes: rmse 0.0118808 max resid 0.1286611
Run 114 stress 0.07520656
Run 115 stress 0.07487575
Run 116 stress 0.0742066
Run 117 stress 0.0741746
Run 118 stress 0.07447713
Run 119 stress 0.07550114
Run 120 stress 0.07531546
Run 121 stress 0.07412308
Run 122 stress 0.07405995
Run 123 stress 0.07536661
Run 124 stress 0.07518678
Run 125 stress 0.0750539
Run 126 stress 0.0742406
Run 127 stress 0.07527225
Run 128 stress 0.07410186
Run 129 stress 0.07486405
Run 130 stress 0.07510652
Run 131 stress 0.07417878
Run 132 stress 0.0743967
Run 133 stress 0.07428291
Run 134 stress 0.0734201
... Procrustes: rmse 0.008337913 max resid 0.09311231
Run 135 stress 0.07494668
Run 136 stress 0.07433737
Run 137 stress 0.07491173
Run 138 stress 0.07506987
Run 139 stress 0.07433888
Run 140 stress 0.07376532
... Procrustes: rmse 0.01132078 max resid 0.1298581
Run 141 stress 0.07426083
Run 142 stress 0.0759338
Run 143 stress 0.07402262
Run 144 stress 0.07419204
Run 145 stress 0.07500366
Run 146 stress 0.07470428
Run 147 stress 0.07505205
Run 148 stress 0.07435951
Run 149 stress 0.07390386
Run 150 stress 0.07399403
Run 151 stress 0.07519878
Run 152 stress 0.07392592
Run 153 stress 0.07406892
Run 154 stress 0.07364193
... Procrustes: rmse 0.008869131 max resid 0.1300762
Run 155 stress 0.07389315
Run 156 stress 0.07494633
Run 157 stress 0.0741264
Run 158 stress 0.07520082
Run 159 stress 0.07526851
Run 160 stress 0.07400545
Run 161 stress 0.07511563
Run 162 stress 0.07518135
Run 163 stress 0.07443978
Run 164 stress 0.07461335
Run 165 stress 0.07505927
Run 166 stress 0.0746276
Run 167 stress 0.07392144
Run 168 stress 0.07497874
Run 169 stress 0.07538232
Run 170 stress 0.07505756
Run 171 stress 0.07522208
Run 172 stress 0.07469917
Run 173 stress 0.07500151
Run 174 stress 0.07389441
Run 175 stress 0.07385946
... Procrustes: rmse 0.01400884 max resid 0.1278604
Run 176 stress 0.07473507
Run 177 stress 0.07414713
Run 178 stress 0.07477712
Run 179 stress 0.07541284
Run 180 stress 0.07393383
Run 181 stress 0.07536193
Run 182 stress 0.07396821
Run 183 stress 0.0744102
Run 184 stress 0.07391726
Run 185 stress 0.0751422
Run 186 stress 0.07405345
Run 187 stress 0.07379454
... Procrustes: rmse 0.01257733 max resid 0.1347662
Run 188 stress 0.07419108
Run 189 stress 0.07504814
Run 190 stress 0.07537465
Run 191 stress 0.07430705
Run 192 stress 0.07392313
Run 193 stress 0.07474243
Run 194 stress 0.0737347
... Procrustes: rmse 0.01018345 max resid 0.1163415
Run 195 stress 0.07393165
Run 196 stress 0.07441157
Run 197 stress 0.07399032
Run 198 stress 0.07382402
... Procrustes: rmse 0.01144517 max resid 0.1288228
Run 199 stress 0.0739376
Run 200 stress 0.07384819
... Procrustes: rmse 0.01054107 max resid 0.1294516
Run 201 stress 0.07496302
Run 202 stress 0.07565388
Run 203 stress 0.07514652
Run 204 stress 0.07498086
Run 205 stress 0.07418292
Run 206 stress 0.07431398
Run 207 stress 0.07387517
Run 208 stress 0.07417817
Run 209 stress 0.07444683
Run 210 stress 0.0745562
Run 211 stress 0.07511673
Run 212 stress 0.07515831
Run 213 stress 0.073975
Run 214 stress 0.07517462
Run 215 stress 0.07571993
Run 216 stress 0.07502441
Run 217 stress 0.07437404
Run 218 stress 0.07460974
Run 219 stress 0.07503461
Run 220 stress 0.07397268
Run 221 stress 0.07411017
Run 222 stress 0.07454477
Run 223 stress 0.07407223
Run 224 stress 0.0739021
Run 225 stress 0.07522509
Run 226 stress 0.0753173
Run 227 stress 0.07384887
... Procrustes: rmse 0.0103963 max resid 0.1140664
Run 228 stress 0.07530704
Run 229 stress 0.07509589
Run 230 stress 0.07522791
Run 231 stress 0.0754115
Run 232 stress 0.07349118
... Procrustes: rmse 0.006143577 max resid 0.09491365
Run 233 stress 0.07404312
Run 234 stress 0.07492119
Run 235 stress 0.07512758
Run 236 stress 0.07402668
Run 237 stress 0.07520954
Run 238 stress 0.07491775
Run 239 stress 0.0755523
Run 240 stress 0.07419907
Run 241 stress 0.07542771
Run 242 stress 0.07423151
Run 243 stress 0.074433
Run 244 stress 0.07507299
Run 245 stress 0.0749457
Run 246 stress 0.07425704
Run 247 stress 0.07542338
Run 248 stress 0.07524835
Run 249 stress 0.07480125
Run 250 stress 0.07424095
Run 251 stress 0.07376875
... Procrustes: rmse 0.009751902 max resid 0.1187984
Run 252 stress 0.07505915
Run 253 stress 0.07394313
Run 254 stress 0.07429427
Run 255 stress 0.07536892
Run 256 stress 0.07511498
Run 257 stress 0.07420732
Run 258 stress 0.07398179
Run 259 stress 0.07429833
Run 260 stress 0.07547262
Run 261 stress 0.07428449
Run 262 stress 0.07429234
Run 263 stress 0.07444704
Run 264 stress 0.07448891
Run 265 stress 0.07415608
Run 266 stress 0.07527221
Run 267 stress 0.07408175
Run 268 stress 0.07509631
Run 269 stress 0.07516579
Run 270 stress 0.07454137
Run 271 stress 0.07541818
Run 272 stress 0.07514797
Run 273 stress 0.07396484
Run 274 stress 0.07521295
Run 275 stress 0.07497417
Run 276 stress 0.07469073
Run 277 stress 0.07373416
... Procrustes: rmse 0.009084333 max resid 0.1202798
Run 278 stress 0.07413188
Run 279 stress 0.07435698
Run 280 stress 0.07533741
Run 281 stress 0.07364552
... Procrustes: rmse 0.01196946 max resid 0.1404214
Run 282 stress 0.07493169
Run 283 stress 0.07467212
Run 284 stress 0.07492725
Run 285 stress 0.07405564
Run 286 stress 0.07433395
Run 287 stress 0.07385769
... Procrustes: rmse 0.01345113 max resid 0.1259407
Run 288 stress 0.07371607
... Procrustes: rmse 0.01199199 max resid 0.1301328
Run 289 stress 0.07391625
Run 290 stress 0.07411003
Run 291 stress 0.07520418
Run 292 stress 0.07394403
Run 293 stress 0.07509057
Run 294 stress 0.0752625
Run 295 stress 0.07539846
Run 296 stress 0.07518962
Run 297 stress 0.07423096
Run 298 stress 0.07393472
Run 299 stress 0.07451568
Run 300 stress 0.07514376
Run 301 stress 0.07434216
Run 302 stress 0.07461602
Run 303 stress 0.07534491
Run 304 stress 0.07519696
Run 305 stress 0.07392339
Run 306 stress 0.07408336
Run 307 stress 0.07406424
Run 308 stress 0.07478821
Run 309 stress 0.07542339
Run 310 stress 0.07371359
... Procrustes: rmse 0.01018135 max resid 0.1307369
Run 311 stress 0.07449296
Run 312 stress 0.07371978
... Procrustes: rmse 0.01009987 max resid 0.1289838
Run 313 stress 0.07386712
Run 314 stress 0.07395433
Run 315 stress 0.07472305
Run 316 stress 0.07481748
Run 317 stress 0.07489642
Run 318 stress 0.07549717
Run 319 stress 0.0750969
Run 320 stress 0.07406464
Run 321 stress 0.07406933
Run 322 stress 0.0744039
Run 323 stress 0.07542437
Run 324 stress 0.07416126
Run 325 stress 0.07527264
Run 326 stress 0.07370265
... Procrustes: rmse 0.01107051 max resid 0.1306372
Run 327 stress 0.0754249
Run 328 stress 0.07431442
Run 329 stress 0.07461884
Run 330 stress 0.07390659
Run 331 stress 0.07462883
Run 332 stress 0.07525901
Run 333 stress 0.07526531
Run 334 stress 0.07424924
Run 335 stress 0.0750202
Run 336 stress 0.07371163
... Procrustes: rmse 0.01242183 max resid 0.1334272
Run 337 stress 0.07387081
Run 338 stress 0.07503239
Run 339 stress 0.07444035
Run 340 stress 0.07369797
... Procrustes: rmse 0.009715762 max resid 0.1364495
Run 341 stress 0.07490495
Run 342 stress 0.07452636
Run 343 stress 0.07514138
Run 344 stress 0.07466526
Run 345 stress 0.07460514
Run 346 stress 0.0751208
Run 347 stress 0.07405165
Run 348 stress 0.07468447
Run 349 stress 0.07545285
Run 350 stress 0.07399098
Run 351 stress 0.07506165
Run 352 stress 0.07545248
Run 353 stress 0.07495447
Run 354 stress 0.07387011
Run 355 stress 0.07458945
Run 356 stress 0.07422489
Run 357 stress 0.07517519
Run 358 stress 0.07422387
Run 359 stress 0.07493579
Run 360 stress 0.07424207
Run 361 stress 0.07466894
Run 362 stress 0.07519356
Run 363 stress 0.07375698
... Procrustes: rmse 0.01079571 max resid 0.1267765
Run 364 stress 0.07502626
Run 365 stress 0.07526004
Run 366 stress 0.07392751
Run 367 stress 0.07507572
Run 368 stress 0.07391212
Run 369 stress 0.07485798
Run 370 stress 0.07406031
Run 371 stress 0.07375768
... Procrustes: rmse 0.01107682 max resid 0.1323648
Run 372 stress 0.0749523
Run 373 stress 0.0739031
Run 374 stress 0.07545146
Run 375 stress 0.07462628
Run 376 stress 0.0750187
Run 377 stress 0.07441731
Run 378 stress 0.0742416
Run 379 stress 0.07516328
Run 380 stress 0.07394378
Run 381 stress 0.07422517
Run 382 stress 0.07438329
Run 383 stress 0.07411848
Run 384 stress 0.07420582
Run 385 stress 0.07501456
Run 386 stress 0.07377241
... Procrustes: rmse 0.01228797 max resid 0.1278678
Run 387 stress 0.07410455
Run 388 stress 0.07519924
Run 389 stress 0.07393327
Run 390 stress 0.07504377
Run 391 stress 0.07561395
Run 392 stress 0.07353481
... Procrustes: rmse 0.009013403 max resid 0.1299546
Run 393 stress 0.07373103
... Procrustes: rmse 0.01231259 max resid 0.1340134
Run 394 stress 0.07513421
Run 395 stress 0.07485496
Run 396 stress 0.07395286
Run 397 stress 0.07398958
Run 398 stress 0.07438072
Run 399 stress 0.07559357
Run 400 stress 0.07475406
Run 401 stress 0.07528713
Run 402 stress 0.07561433
Run 403 stress 0.07474959
Run 404 stress 0.07460552
Run 405 stress 0.07473478
Run 406 stress 0.07394682
Run 407 stress 0.07478529
Run 408 stress 0.07514026
Run 409 stress 0.07419717
Run 410 stress 0.07393309
Run 411 stress 0.07430879
Run 412 stress 0.07397659
Run 413 stress 0.07358595
... Procrustes: rmse 0.01306045 max resid 0.1258173
Run 414 stress 0.074179
Run 415 stress 0.07480772
Run 416 stress 0.07440817
Run 417 stress 0.07555891
Run 418 stress 0.07403137
Run 419 stress 0.07467252
Run 420 stress 0.07481558
Run 421 stress 0.0744667
Run 422 stress 0.07475577
Run 423 stress 0.07537888
Run 424 stress 0.07512247
Run 425 stress 0.07419997
Run 426 stress 0.07348721
... Procrustes: rmse 0.01223881 max resid 0.1276855
Run 427 stress 0.07408208
Run 428 stress 0.07399277
Run 429 stress 0.07443146
Run 430 stress 0.07525452
Run 431 stress 0.07369938
... Procrustes: rmse 0.01272754 max resid 0.1268051
Run 432 stress 0.07598447
Run 433 stress 0.07557919
Run 434 stress 0.07390075
Run 435 stress 0.07530529
Run 436 stress 0.07499659
Run 437 stress 0.07528059
Run 438 stress 0.07490328
Run 439 stress 0.07362399
... Procrustes: rmse 0.01292832 max resid 0.1259931
Run 440 stress 0.07394942
Run 441 stress 0.07407381
Run 442 stress 0.07420546
Run 443 stress 0.07575533
Run 444 stress 0.07468921
Run 445 stress 0.07482684
Run 446 stress 0.07485404
Run 447 stress 0.07434552
Run 448 stress 0.07402306
Run 449 stress 0.07393731
Run 450 stress 0.07518266
Run 451 stress 0.07397705
Run 452 stress 0.07401406
Run 453 stress 0.07504258
Run 454 stress 0.0749351
Run 455 stress 0.07414586
Run 456 stress 0.07565358
Run 457 stress 0.07438571
Run 458 stress 0.07412992
Run 459 stress 0.07519675
Run 460 stress 0.0740214
Run 461 stress 0.07438382
Run 462 stress 0.0755578
Run 463 stress 0.07441952
Run 464 stress 0.07446079
Run 465 stress 0.07414533
Run 466 stress 0.07516437
Run 467 stress 0.07409466
Run 468 stress 0.07418167
Run 469 stress 0.0740646
Run 470 stress 0.07522238
Run 471 stress 0.07455339
Run 472 stress 0.07601397
Run 473 stress 0.07402227
Run 474 stress 0.07531606
Run 475 stress 0.0744596
Run 476 stress 0.07509474
Run 477 stress 0.0748124
Run 478 stress 0.07530191
Run 479 stress 0.07522585
Run 480 stress 0.07493452
Run 481 stress 0.07450819
Run 482 stress 0.07453053
Run 483 stress 0.07528691
Run 484 stress 0.07510803
Run 485 stress 0.07379753
... Procrustes: rmse 0.009806388 max resid 0.1298822
Run 486 stress 0.07544587
Run 487 stress 0.07397591
Run 488 stress 0.0748057
Run 489 stress 0.07371733
... Procrustes: rmse 0.01116551 max resid 0.1282052
Run 490 stress 0.07435777
Run 491 stress 0.07440119
Run 492 stress 0.07528999
Run 493 stress 0.07468488
Run 494 stress 0.07421266
Run 495 stress 0.07458019
Run 496 stress 0.07352401
... Procrustes: rmse 0.006242808 max resid 0.08835656
Run 497 stress 0.07385832
... Procrustes: rmse 0.01014896 max resid 0.1302035
Run 498 stress 0.07417513
Run 499 stress 0.07425166
Run 500 stress 0.075162
*** No convergence -- monoMDS stopping criteria:
500: scale factor of the gradient < sfgrmin
Still no convergence.
That looks pretty damn good though.
Also not bad.
Let’s go forward with that for the time being.
Plotting the NMDS
# Plotting points in ordination space
plot(NMS, "sites") # Produces distance
orditorp(NMS, "sites") # Gives points labels
BUT, I’m actually less interested in this and more interested in using it with the sample replicates.
Here is where the uncollapsed df might be more helpful… let me take a quick look to see what Pat did here (and if what I did for Cornell might be useful)